library(R2jags)
library(foreign)

setwd("c:/users/james/desktop/dropbox/transparency_and_democracy/transparencyindex/PAReplicationMaterials/regressions")

load("TransparencyABPrep.RData")

attach.jags(results)

k<-0
transparency.means<-c(1:2120)

for (i in 1:106){
  for (j in 1:20){
  k<-k+1
  transparency.means[k] <- mean(transparency[,i,j])
  }
}

HRVmargeffectCorrup.sim <- array(NA, c(2000,100))
PapermargeffectCorrup.sim<- array(NA, c(2000,100))
HRVmargeffectLawOrd.sim <- array(NA, c(2000,100))
PapermargeffectLawOrd.sim<- array(NA, c(2000,100))
HRVmargeffectBurQual.sim <- array(NA, c(2000,100))
PapermargeffectBurQual.sim<- array(NA, c(2000,100))

polity.sim<- -10

for (i in 1:100){
  HRVmargeffectCorrup.sim[, i]  <-  gamma.corrup[,4]*sd(transparency.means) + gamma.corrup[,5]*sd(transparency.means)*polity.sim
  PapermargeffectCorrup.sim[,i] <- gamma.corrup[,6]*sd(newscirc) + gamma.corrup[,7]*sd(newscirc)*polity.sim
  HRVmargeffectLawOrd.sim[, i]  <-  gamma.laword[,4]*sd(transparency.means) + gamma.laword[,5]*sd(transparency.means)*polity.sim
  PapermargeffectLawOrd.sim[,i] <- gamma.laword[,6]*sd(newscirc) + gamma.laword[,7]*sd(newscirc)*polity.sim
  HRVmargeffectBurQual.sim[, i]  <-  gamma.burqual[,4]*sd(transparency.means) + gamma.burqual[,5]*sd(transparency.means)*polity.sim
  PapermargeffectBurQual.sim[,i] <- gamma.burqual[,6]*sd(newscirc) + gamma.burqual[,7]*sd(newscirc)*polity.sim
  polity.sim<- polity.sim + .2
  }
  
HRVcorrupEst <- c(1:100)

HRVcorrupLB <- c(1:100)
HRVcorrupUB <- c(1:100)

HRVcorrupHPD <- cbind(HRVcorrupLB, HRVcorrupUB)

PapercorrupEst <- c(1:100)

PapercorrupLB <- c(1:100)
PapercorrupUB <- c(1:100)

PapercorrupHPD <- cbind(PapercorrupLB, PapercorrupUB)

HRVlawordEst <- c(1:100)

HRVlawordLB <- c(1:100)
HRVlawordUB <- c(1:100)

HRVlawordHPD <- cbind(HRVlawordLB, HRVlawordUB)

PaperlawordEst <- c(1:100)

PaperlawordLB <- c(1:100)
PaperlawordUB <- c(1:100)

PaperlawordHPD <- cbind(PaperlawordLB, PaperlawordUB)

HRVburqualEst <- c(1:100)

HRVburqualLB <- c(1:100)
HRVburqualUB <- c(1:100)

HRVburqualHPD <- cbind(HRVburqualLB, HRVburqualUB)

PaperburqualEst <- c(1:100)

PaperburqualLB <- c(1:100)
PaperburqualUB <- c(1:100)

PaperburqualHPD <- cbind(PaperburqualLB, PaperburqualUB)

for (i in 1:100){
  HRVcorrupEst[i] <- mean(HRVmargeffectCorrup.sim[,i])
  HRVlawordEst[i] <- mean(HRVmargeffectLawOrd.sim[,i])
  HRVburqualEst[i] <- mean(HRVmargeffectBurQual.sim[,i])

  PapercorrupEst[i] <- mean(PapermargeffectCorrup.sim[,i])
  PaperlawordEst[i] <- mean(PapermargeffectLawOrd.sim[,i])
  PaperburqualEst[i] <- mean(PapermargeffectBurQual.sim[,i])
  
  HRVcorrupHPD[i,] <- HPDinterval(as.mcmc(as.matrix(HRVmargeffectCorrup.sim[,i])))
  HRVlawordHPD[i,] <- HPDinterval(as.mcmc(as.matrix(HRVmargeffectLawOrd.sim[,i])))
  HRVburqualHPD[i,] <- HPDinterval(as.mcmc(as.matrix(HRVmargeffectBurQual.sim[,i])))
  
  PapercorrupHPD[i,] <- HPDinterval(as.mcmc(as.matrix(PapermargeffectCorrup.sim[,i])))
  PaperlawordHPD[i,] <- HPDinterval(as.mcmc(as.matrix(PapermargeffectLawOrd.sim[,i])))
  PaperburqualHPD[i,] <- HPDinterval(as.mcmc(as.matrix(PapermargeffectBurQual.sim[,i])))
  }
  
polity.sim<-c(1:100)
polity.sim<-(polity.sim-51)/5

pdf('HRVCorruptEffect.pdf')
plot(x=c(), y=c(), xlim=c(-10, 10), ylim=c(-1,1), xlab="Polity2", ylab="Marginal Corruption Effect", main="HRV Scores and Corruption")
lines(y=HRVcorrupEst, x=polity.sim)
lines(y=HRVcorrupHPD[,1], x=polity.sim, lty=2)
lines(y=HRVcorrupHPD[,2], x=polity.sim, lty=2)
abline(h=0, lty=3)
dev.off()


pdf('HRVlawordEffect.pdf')
plot(x=c(), y=c(), xlim=c(-10, 10), ylim=c(-1,1), xlab="Polity2", ylab="Marginal Corruption Effect", main="HRV Scores and Law and Order")
lines(y=HRVlawordEst, x=polity.sim)
lines(y=HRVlawordHPD[,1], x=polity.sim, lty=2)
lines(y=HRVlawordHPD[,2], x=polity.sim, lty=2)
abline(h=0, lty=3)
dev.off()

pdf('HRVburqualEffect.pdf')
plot(x=c(), y=c(), xlim=c(-10, 10), ylim=c(-1,1), xlab="Polity2", ylab="Marginal Corruption Effect", main="HRV Scores and Bureaucratic Quality")
lines(y=HRVburqualEst, x=polity.sim)
lines(y=HRVburqualHPD[,1], x=polity.sim, lty=2)
lines(y=HRVburqualHPD[,2], x=polity.sim, lty=2)
abline(h=0, lty=3)
dev.off()

pdf('PaperCorruptEffect.pdf')
plot(x=c(), y=c(), xlim=c(-10, 10), ylim=c(-1,1), xlab="Polity2", ylab="Marginal Corruption Effect", main="News Circulation and Corruption")
lines(y=PapercorrupEst, x=polity.sim)
lines(y=PapercorrupHPD[,1], x=polity.sim, lty=2)
lines(y=PapercorrupHPD[,2], x=polity.sim, lty=2)
abline(h=0, lty=3)
dev.off()


pdf('PaperlawordEffect.pdf')
plot(x=c(), y=c(), xlim=c(-10, 10), ylim=c(-1,1), xlab="Polity2", ylab="Marginal Corruption Effect", main="News Circulation Scores and Law and Order")
lines(y=PaperlawordEst, x=polity.sim)
lines(y=PaperlawordHPD[,1], x=polity.sim, lty=2)
lines(y=PaperlawordHPD[,2], x=polity.sim, lty=2)
abline(h=0, lty=3)
dev.off()

pdf('PaperburqualEffect.pdf')
plot(x=c(), y=c(), xlim=c(-10, 10), ylim=c(-1,1), xlab="Polity2", ylab="Marginal Corruption Effect", main="News Circulation and Bureaucratic Quality")
lines(y=PaperburqualEst, x=polity.sim)
lines(y=PaperburqualHPD[,1], x=polity.sim, lty=2)
lines(y=PaperburqualHPD[,2], x=polity.sim, lty=2)
abline(h=0, lty=3)
dev.off()
